| Variable | Definition |
|---|---|
| K | Fulton’s Condition Factor (<1=good, >1=bad) |
| W | weight of host fish (in grams) |
| SL | standard length of fish (in millimeters) |
\[ K = W/SL^3 *100 \]
library(ggeffects)
library(tidyverse)
library(ggpubr)
library(caret)
library(tibble)
library(readr)
library(dbplyr)
library(lubridate)
library(ggplot2)
library(glmmTMB)
Pvig.life.stage <- P.vigilax.data %>%
mutate(adult = rowSums(select(.,
"acanth_bigb",
"acanth_bkr",
"nem_cona",
"nem_cap",
"nem_dich",
"nem_larv",
"nem_trbk",
"nem_unk",
"cest_godz",
"cest_comp",
"cest_ntg",
"cope_cl",
"cope_grab",
"cest_pea",
"cest_l",
"cest_vit",
"cest_unk",
"mono_gyro",
"mono_dact",
"mono_lg",
"trem_rnda",
"trem_unk"))) %>%
mutate(larval = rowSums(select(.,
"meta_hs",
"meta_unk",
"myx_go",
"myx_geom",
"myxo_bc",
"myx_thel",
"myxo_sbad",
"myxo_unk",
"myxo_up",
"trem_met",
"trem_sss",
"trem_buc",
"trem_metaf",
"trem_metags",
"trem_ms",
"trem_nea"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
janitor::clean_names()
Ipun.life.stage <- I.punctatus.data %>%
mutate(larval = rowSums(select(.,
"nem_cyst",
"trem_meta_gg",
"trem_meta_unk",
"trem_lg",
"nem_cyst",
"myx_tail",
"myx_geom"))) %>%
mutate(adult = rowSums(select(.,
"acan_ostr",
"cest_comp",
"cest_meg",
"leech_kut",
"nem_nmts",
"mono_ip",
"mono_unk",
"nem_phar",
"nem_pty",
"nem_sp",
"nem_taf",
"nem_unk",
"nmorph",
"trem_2",
"trem_allo",
"trem_ble",
"trem_crep",
"trem_f",
"trem_lip",
"trem_megict",
"trem_smile",
"trem_unk"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
janitor::clean_names()
Nath.life.stage <- N.atherinoides.data %>%
mutate(larval = rowSums(select(.,
"nem_cyst",
"trem_larv",
"trem_meta",
"myx_sp",
"myx_round",
"trem_cyst"))) %>%#
mutate(adult = rowSums(select(.,
"cest_comp",
"cest_hyd",
"cest_lvs",
"cest_tri",
"cope_a",
"mono_all",
"mono_unk",
"nem_bran",
"nem_cap",
"nem_spky",
"nem_tfs",
"nem_twas",
"nem_unkn",
"trem_gb",
"trem_l"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
janitor::clean_names()
####
Pvig.life.stage <- P.vigilax.data %>%
mutate(adult = rowSums(select(.,
"acanth_bigb",
"acanth_bkr",
"nem_cona",
"nem_cap",
"nem_dich",
"nem_larv",
"nem_trbk",
"nem_unk",
"cest_godz",
"cest_comp",
"cest_ntg",
"cope_cl",
"cope_grab",
"cest_pea",
"cest_l",
"cest_vit",
"cest_unk",
"mono_gyro",
"mono_dact",
"mono_lg",
"trem_rnda",
"trem_unk"))) %>%
mutate(larval = rowSums(select(.,
"meta_hs",
"meta_unk",
"myx_go",
"myx_geom",
"myxo_bc",
"myx_thel",
"myxo_sbad",
"myxo_unk",
"myxo_up",
"trem_met",
"trem_sss",
"trem_buc",
"trem_metaf",
"trem_metags",
"trem_ms",
"trem_nea"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
select(-matches("^(myx_|myxo_|acanth_|mono_|cope_)")) %>%
janitor::clean_names()
####
Ipun.life.stage <- I.punctatus.data %>%
mutate(larval = rowSums(select(.,
"nem_cyst",
"trem_meta_gg",
"trem_meta_unk",
"trem_lg",
"nem_cyst",
"myx_tail",
"myx_geom"))) %>%
mutate(adult = rowSums(select(.,
"acan_ostr",
"cest_comp",
"cest_meg",
"leech_kut",
"nem_nmts",
"mono_ip",
"mono_unk",
"nem_phar",
"nem_pty",
"nem_sp",
"nem_taf",
"nem_unk",
"nmorph",
"trem_2",
"trem_allo",
"trem_ble",
"trem_crep",
"trem_f",
"trem_lip",
"trem_megict",
"trem_smile",
"trem_unk"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
select(-matches("^(myx_|acan_|mono_|leech_|nmorph)")) %>%
janitor::clean_names()
####
Nath.life.stage <- N.atherinoides.data %>%
mutate(larval = rowSums(select(.,
"nem_cyst",
"trem_larv",
"trem_meta",
"myx_sp",
"myx_round",
"trem_cyst"))) %>%#
mutate(adult = rowSums(select(.,
"cest_comp",
"cest_hyd",
"cest_lvs",
"cest_tri",
"cope_a",
"mono_all",
"mono_unk",
"nem_bran",
"nem_cap",
"nem_spky",
"nem_tfs",
"nem_twas",
"nem_unkn",
"trem_gb",
"trem_l"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
select(-matches("^(myx_|mono_|cope_)")) %>%
janitor::clean_names()
Nath.life.stage <- Nath.life.stage %>%
mutate(body.index.ath = ((Nath.life.stage$weight_mg/((Nath.life.stage$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
Pvig.life.stage <- Pvig.life.stage %>%
mutate(body.index.vig = ((Pvig.life.stage$weight_mg/((Pvig.life.stage$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
Ipun.life.stage <- Ipun.life.stage %>%
mutate(body.index.pun = ((Ipun.life.stage$weight_mg/((Ipun.life.stage$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
na_ip_pv_combined <- bind_rows(Nath.life.stage, Pvig.life.stage, Ipun.life.stage, .id = "dataset")
na_ip_pv_pp <- preProcess(na_ip_pv_combined, method = c("center", "scale"))
na_ip_pv_scaled <- predict(na_ip_pv_pp, na_ip_pv_combined)
notath_scaled <- na_ip_pv_scaled %>% filter(dataset == "1") %>%
mutate(body.index.notath = ((Nath.life.stage$weight_mg/((Nath.life.stage$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
pimvig_scaled <- na_ip_pv_scaled %>% filter(dataset == "2") %>%
mutate(body.index.pimvig = ((Pvig.life.stage$weight_mg/((Pvig.life.stage$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
ictpun_scaled <- na_ip_pv_scaled %>% filter(dataset == "3") %>%
mutate(body.index.ictpun = ((Ipun.life.stage$weight_mg/((Ipun.life.stage$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
# plot condition factor vs binary code for life stage (aka "stage-specific parasite stuff")
NA.bc.a <- ggplot(Nath.life.stage,
aes(x = adult, y = body_index_ath, color = sex)) + xlab("Adult") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
NA.bc.l <- ggplot(Nath.life.stage,
aes(x = larval, y = body_index_ath, color = sex)) + xlab("Larval") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
# plot condition factor vs binary code for life stage (aka "stage-specific parasite stuff")
IP.bc.a <- ggplot(Ipun.life.stage,
aes(x = adult, y = body_index_pun, color = sex)) + xlab("Adult") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
IP.bc.l <- ggplot(Ipun.life.stage,
aes(x = larval, y = body_index_pun, color = sex)) + xlab("Larval") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
# plot condition factor vs binary code for life stage (aka "stage-specific parasite stuff")
PV.bc.a <- ggplot(Pvig.life.stage,
aes(x = adult, y = body_index_vig, color = sex)) + xlab("Adult") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
PV.bc.l <- ggplot(Pvig.life.stage,
aes(x = larval, y = body_index_vig, color = sex)) + xlab("Larval") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
ggarrange(IP.bc.a, IP.bc.l, NA.bc.a, NA.bc.l, PV.bc.a, PV.bc.a, common.legend = TRUE, legend = "bottom", labels = c("A", "B", "C", "D", "E", "F"))
# more variability than the others
## double check sex, weights
## should we exclude these outliers
# plot condition factor vs binary code for life stage (aka "stage-specific parasite stuff")
NA.bc.a <- ggplot(notath_scaled,
aes(x = adult, y = body_index_ath, color = sex)) + xlab("Adult") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
NA.bc.l <- ggplot(notath_scaled,
aes(x = larval, y = body_index_ath, color = sex)) + xlab("Larval") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
# plot condition factor vs binary code for life stage (aka "stage-specific parasite stuff")
IP.bc.a <- ggplot(ictpun_scaled,
aes(x = adult, y = body_index_pun, color = sex)) + xlab("Adult") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
IP.bc.l <- ggplot(ictpun_scaled,
aes(x = larval, y = body_index_pun, color = sex)) + xlab("Larval") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
# plot condition factor vs binary code for life stage (aka "stage-specific parasite stuff")
PV.bc.a <- ggplot(pimvig_scaled,
aes(x = adult, y = body_index_vig, color = sex)) + xlab("Adult") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
PV.bc.l <- ggplot(pimvig_scaled,
aes(x = larval, y = body_index_vig, color = sex)) + xlab("Larval") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
ggarrange(IP.bc.a, IP.bc.l, NA.bc.a, NA.bc.l, PV.bc.a, PV.bc.a, common.legend = TRUE, legend = "bottom", labels = c("A", "B", "C", "D", "E", "F"))
# more variability than the others
## double check sex, weights
## should we exclude these outliers
We want to run a general linear model to test the research questions: 1. Is body condition associated with increased total parasite burden? - is there a linear association between increasing body condition and increasing parasite burden? - does FCF increase as parasite burden increases? 2. Is body condition associated with infection with adult or larval stages of parasites? - is there a linear association between increasing body condition and increasing abundance of larval or adult parasites in the host? - does FCF increase as adult or larval parasite presence increases? - variables: - predictor = parasite_count, life_stage –> discrete - response = FCF –> continuous - discrete predictor = ANOVA - multiple/ continuous and discrete predictors = ANCOVA?
\[ K_i\sim\alpha+\beta_1count+\epsilon_i \]
\[ K_i\sim\alpha+\beta_1stage+\beta_2count+\epsilon_i \] ### Run linear model to test above questions
# Question 1 (lm)
summary(FCF_ath_lm <- lm(body_index_ath ~ psite_count, data = Nath.life.stage))
##
## Call:
## lm(formula = body_index_ath ~ psite_count, data = Nath.life.stage)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.8708 -0.2325 -0.1553 -0.0608 7.2784
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.088027 0.085561 12.716 <2e-16 ***
## psite_count 0.001006 0.003980 0.253 0.801
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9823 on 169 degrees of freedom
## (37 observations deleted due to missingness)
## Multiple R-squared: 0.0003782, Adjusted R-squared: -0.005537
## F-statistic: 0.06395 on 1 and 169 DF, p-value: 0.8007
plot(FCF_ath_lm)
summary(FCF_ath_lm <- lm(body_index_pun ~ psite_count, data = Ipun.life.stage))
##
## Call:
## lm(formula = body_index_pun ~ psite_count, data = Ipun.life.stage)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.8625 -0.5370 -0.3616 -0.1666 13.7487
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.860079 0.261188 7.122 6.59e-10 ***
## psite_count -0.002641 0.005931 -0.445 0.657
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.14 on 72 degrees of freedom
## (15 observations deleted due to missingness)
## Multiple R-squared: 0.002746, Adjusted R-squared: -0.01111
## F-statistic: 0.1982 on 1 and 72 DF, p-value: 0.6575
plot(FCF_ath_lm)
summary(FCF_ath_lm <- lm(body_index_vig ~ psite_count, data = Pvig.life.stage))
##
## Call:
## lm(formula = body_index_vig ~ psite_count, data = Pvig.life.stage)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.8966 -0.4336 -0.2843 -0.0619 12.7554
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.997946 0.305819 6.533 3.51e-08 ***
## psite_count -0.003891 0.008464 -0.460 0.648
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.874 on 49 degrees of freedom
## (146 observations deleted due to missingness)
## Multiple R-squared: 0.004294, Adjusted R-squared: -0.01603
## F-statistic: 0.2113 on 1 and 49 DF, p-value: 0.6478
plot(FCF_ath_lm)
# Question 1 (ANOVA)
summary(FCF_ath_anova <- aov(body_index_ath ~ psite_count, data = Nath.life.stage))
## Df Sum Sq Mean Sq F value Pr(>F)
## psite_count 1 0.06 0.0617 0.064 0.801
## Residuals 169 163.07 0.9649
## 37 observations deleted due to missingness
plot(FCF_ath_anova)
summary(FCF_ath_anova <- aov(body_index_pun ~ psite_count, data = Ipun.life.stage))
## Df Sum Sq Mean Sq F value Pr(>F)
## psite_count 1 0.9 0.907 0.198 0.657
## Residuals 72 329.6 4.578
## 15 observations deleted due to missingness
plot(FCF_ath_anova)
summary(FCF_ath_anova <- aov(body_index_vig ~ psite_count, data = Pvig.life.stage))
## Df Sum Sq Mean Sq F value Pr(>F)
## psite_count 1 0.74 0.742 0.211 0.648
## Residuals 49 172.13 3.513
## 146 observations deleted due to missingness
plot(FCF_ath_anova)
# Question 2 (lm)
summary(FCF_ath_al_lm <- lm(body_index_ath ~ larval * adult * psite_count, data = Nath.life.stage))
##
## Call:
## lm(formula = body_index_ath ~ larval * adult * psite_count, data = Nath.life.stage)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.9419 -0.2419 -0.1293 -0.0078 7.0639
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.017e+00 1.095e-01 9.290 <2e-16 ***
## larval 4.827e-03 1.537e-02 0.314 0.754
## adult 1.220e-02 2.597e-02 0.470 0.639
## psite_count NA NA NA NA
## larval:adult 5.821e-03 3.952e-03 1.473 0.143
## larval:psite_count -3.455e-05 9.992e-05 -0.346 0.730
## adult:psite_count -5.134e-04 6.986e-04 -0.735 0.463
## larval:adult:psite_count -6.880e-05 5.109e-05 -1.347 0.180
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9831 on 164 degrees of freedom
## (37 observations deleted due to missingness)
## Multiple R-squared: 0.02828, Adjusted R-squared: -0.007275
## F-statistic: 0.7954 on 6 and 164 DF, p-value: 0.5748
plot(FCF_ath_al_lm)
summary(FCF_ath_al_lm <- lm(body_index_pun ~ larval * adult * psite_count, data = Ipun.life.stage))
##
## Call:
## lm(formula = body_index_pun ~ larval * adult * psite_count, data = Ipun.life.stage)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.9417 -0.5111 -0.3990 -0.1282 13.6316
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.9771812 0.3748861 5.274 1.54e-06 ***
## larval -0.2912635 0.5894503 -0.494 0.623
## adult -0.0220916 0.0757006 -0.292 0.771
## psite_count NA NA NA NA
## larval:adult 0.0190266 0.0422251 0.451 0.654
## larval:psite_count 0.0040358 0.0079439 0.508 0.613
## adult:psite_count 0.0002453 0.0016887 0.145 0.885
## larval:adult:psite_count -0.0004097 0.0008118 -0.505 0.615
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.209 on 67 degrees of freedom
## (15 observations deleted due to missingness)
## Multiple R-squared: 0.01098, Adjusted R-squared: -0.07758
## F-statistic: 0.124 on 6 and 67 DF, p-value: 0.9931
plot(FCF_ath_al_lm)
summary(FCF_ath_al_lm <- lm(body_index_vig ~ larval * adult * psite_count, data = Pvig.life.stage))
##
## Call:
## lm(formula = body_index_vig ~ larval * adult * psite_count, data = Pvig.life.stage)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.9284 -0.4533 -0.2908 -0.0267 12.6762
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.104e+00 4.792e-01 4.390 7.01e-05 ***
## larval -2.445e-02 5.497e-02 -0.445 0.659
## adult 4.899e-02 3.044e-01 0.161 0.873
## psite_count NA NA NA NA
## larval:adult 1.423e-02 4.344e-02 0.327 0.745
## larval:psite_count 1.495e-04 3.699e-04 0.404 0.688
## adult:psite_count -1.132e-02 3.694e-02 -0.306 0.761
## larval:adult:psite_count -2.054e-05 9.178e-05 -0.224 0.824
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.972 on 44 degrees of freedom
## (146 observations deleted due to missingness)
## Multiple R-squared: 0.01026, Adjusted R-squared: -0.1247
## F-statistic: 0.07602 on 6 and 44 DF, p-value: 0.9981
plot(FCF_ath_al_lm)
# Question 2 (ANOVA)
summary(FCF_ath_al_anova <- aov(body_index_ath ~ larval * adult * psite_count, data = Nath.life.stage))
## Df Sum Sq Mean Sq F value Pr(>F)
## larval 1 0.01 0.0079 0.008 0.928
## adult 1 0.15 0.1524 0.158 0.692
## larval:adult 1 0.01 0.0100 0.010 0.919
## larval:psite_count 1 0.80 0.8010 0.829 0.364
## adult:psite_count 1 1.89 1.8884 1.954 0.164
## larval:adult:psite_count 1 1.75 1.7530 1.814 0.180
## Residuals 164 158.52 0.9666
## 37 observations deleted due to missingness
plot(FCF_ath_al_anova)
summary(FCF_ath_al_anova <- aov(body_index_pun ~ larval * adult * psite_count, data = Ipun.life.stage))
## Df Sum Sq Mean Sq F value Pr(>F)
## larval 1 0.5 0.532 0.109 0.742
## adult 1 1.1 1.054 0.216 0.644
## larval:adult 1 0.5 0.470 0.096 0.757
## larval:psite_count 1 0.2 0.245 0.050 0.823
## adult:psite_count 1 0.1 0.086 0.018 0.895
## larval:adult:psite_count 1 1.2 1.243 0.255 0.615
## Residuals 67 326.9 4.879
## 15 observations deleted due to missingness
plot(FCF_ath_al_anova)
summary(FCF_ath_al_anova <- aov(body_index_vig ~ larval * adult * psite_count, data = Pvig.life.stage))
## Df Sum Sq Mean Sq F value Pr(>F)
## larval 1 0.67 0.669 0.172 0.680
## adult 1 0.31 0.305 0.079 0.781
## larval:adult 1 0.01 0.007 0.002 0.967
## larval:psite_count 1 0.37 0.375 0.096 0.758
## adult:psite_count 1 0.22 0.223 0.057 0.812
## larval:adult:psite_count 1 0.19 0.195 0.050 0.824
## Residuals 44 171.10 3.889
## 146 observations deleted due to missingness
plot(FCF_ath_al_anova)